---
title: "gender_protesting_2023_analysis"
output: pdf_document
date: "2023-10-01"
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning=FALSE,message=FALSE)
setwd("/Users/martinnaunov/Desktop/Gender Protesting Data")
```

```{r, include=FALSE}
rm(list=ls())
```

```{r}
setwd("/Users/martinnaunov/Desktop/Gender Protesting Data")
gp_1 <- read.csv("gp_c.csv")

```


#MAIN ANALYSIS

```{r}
gp_1$condition_f <- as.factor(gp_1$condition)

```

```{r}
v_1_model <- aov(violence_1 ~ condition_f, data = gp_1)

v1hsd <- TukeyHSD(v_1_model)
```


```{r}

v_2_model <- aov(violence_2 ~ condition_f, data = gp_1)

v2hsd <- TukeyHSD(v_2_model)
```

```{r}

t_1_model <- aov(toler_1_index ~ condition_f, data = gp_1)

t1hsd <- TukeyHSD(t_1_model)
```

```{r}

t_2_model <- aov(toler_2_index ~ condition_f, data = gp_1)

t2hsd <- TukeyHSD(t_2_model)
```

```{r}

m_model <- aov(mor_1 ~ condition_f, data = gp_1)
mhsd <- TukeyHSD(m_model)

```


#Figure 3 (Differences of Means for Repression Toleration)

```{r}
library(broom)
t1hsd_df<- tidy(t1hsd)
t2hsd_df<- tidy(t2hsd)
```

```{r}
library(dplyr)

t1hsd_df_2 <- t1hsd_df[c("contrast", "estimate","conf.low","conf.high", "adj.p.value")]
t2hsd_df_2 <- t2hsd_df[c("contrast", "estimate","conf.low","conf.high", "adj.p.value")]

t1hsd_df_2$dv <- "1"
t2hsd_df_2$dv <- "2"

hsd_df_t1t2 <- rbind(t1hsd_df_2, t2hsd_df_2)

hsd_df_t1t2_t <-hsd_df_t1t2[!(hsd_df_t1t2$contrast=="3-2" | hsd_df_t1t2$contrast=="4-2"),]

i1 <- hsd_df_t1t2_t$contrast=="4-3"
hsd_df_t1t2_t$estimate[i1] <- hsd_df_t1t2_t$estimate[i1] * (-1)
hsd_df_t1t2_t$conf.low[i1] <- hsd_df_t1t2_t$conf.low[i1] * (-1)
hsd_df_t1t2_t$conf.high[i1] <- hsd_df_t1t2_t$conf.high[i1] * (-1)

```


```{r}

plot_repression_f3 <- hsd_df_t1t2_t
```


```{r}
library(ggplot2)
library(dplyr)

plot_repression_f3$names <- plot_repression_f3$dv

plot_repression_f3$names<-dplyr::recode(plot_repression_f3$names, 
                          "1"="Preemptive \nRepression \nToleration",
                          "2"="Reactive \nRepression \nToleration")

plot_repression_f3$namesf<-as.factor(plot_repression_f3$names)

plot_repression_f3$namesfd <- factor(plot_repression_f3$namesf, levels = c("Reactive \nRepression \nToleration","Preemptive \nRepression \nToleration"))

plot_repression_f3$group <- plot_repression_f3$contrast

plot_repression_f3$group <-dplyr::recode(plot_repression_f3$group, "2-1"="Generic Women \nvs. Men", "3-1"="Patriarchy-Compliant \nWomen vs. Men", "4-1"="Patriarchy-Defiant \nWomen vs. Men", "4-3"="Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women")

plot_repression_f3$group_f<-as.factor(plot_repression_f3$group)

plot_repression_f3$group_f <- factor(plot_repression_f3$group_f, levels = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "Patriarchy-Defiant \nWomen vs. Men", "Patriarchy-Compliant \nWomen vs. Men", "Generic Women \nvs. Men"))

plot_t1t2_f3<-ggplot(plot_repression_f3, aes(x = estimate, y = namesfd, group=group_f, color=group_f))+ geom_point(aes(shape=group_f, color=group_f), size=3, position = position_dodge(width = 0.50))+ geom_errorbar(aes(xmin=conf.low, xmax=conf.high, group=group_f), size=1, width=0.3,position = position_dodge(width = 0.50)) + geom_vline(xintercept = 0,lty=2) + labs(title = " ", subtitle = "") + theme_bw()+theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5), axis.text.y = element_text(color = "grey20", size = 12, face = "bold"),legend.text=element_text(size=12, face = "bold")) + xlab("") + ylab("")


pl_f3 <- plot_t1t2_f3 + labs(color = " ", shape=" ") + scale_colour_grey(start = 0.15, end = 0.85)


p_f3 <- pl_f3 + theme(
  legend.position = "bottom", 
  legend.justification = "right", # Move the legend to the left
  legend.margin = margin(r = 7, unit = "mm") # Adjust this value to control the distance from the left side
)

p_f3

p_f3_c <- p_f3 +
  scale_shape_manual(values = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women" = 15, 
                                "Patriarchy-Defiant \nWomen vs. Men" = 17,  
                                "Patriarchy-Compliant \nWomen vs. Men" = 4,  
                                "Generic Women \nvs. Men" = 16)) 



# Reverse the order of the legend items
plot_t1t2_f3 <- p_f3_c + guides(color = guide_legend(reverse = TRUE),
                                     shape = guide_legend(reverse = TRUE))


plot_t1t2_f3


ggsave("GP_Figure3.tiff", plot_t1t2_f3, width = 9, height = 6, dpi = 1000,
       units = "in")


```



#Figure 4 (Differences of Means for Perception of Protest Violence)

```{r}
library(broom)
v1hsd_df<- tidy(v1hsd)
v2hsd_df<- tidy(v2hsd)
```

```{r}
v1hsd_df_2 <- v1hsd_df[c("contrast", "estimate","conf.low","conf.high", "adj.p.value")]
v2hsd_df_2 <- v2hsd_df[c("contrast", "estimate","conf.low","conf.high", "adj.p.value")]

v1hsd_df_2$dv <- "1"
v2hsd_df_2$dv <- "2"

hsd_df_v1v2 <- rbind(v1hsd_df_2, v2hsd_df_2)

hsd_df_v1v2_t <-hsd_df_v1v2[!(hsd_df_v1v2$contrast=="3-2" | hsd_df_v1v2$contrast=="4-2"),]

i1 <- hsd_df_v1v2_t$contrast=="4-3"
hsd_df_v1v2_t$estimate[i1] <- hsd_df_v1v2_t$estimate[i1] * (-1)
hsd_df_v1v2_t$conf.low[i1] <- hsd_df_v1v2_t$conf.low[i1] * (-1)
hsd_df_v1v2_t$conf.high[i1] <- hsd_df_v1v2_t$conf.high[i1] * (-1)
```

```{r}
plot_violence_f4 <- hsd_df_v1v2_t
```

```{r}
library(ggplot2)

#DV
plot_violence_f4$names <- plot_violence_f4$dv

plot_violence_f4$names<-dplyr::recode(plot_violence_f4$names, 
                          "1"="Preemptive \nViolence \nPerception",
                          "2"="Believing Authorities \nthat Protesters \nwere Violent")

plot_violence_f4$namesf<-as.factor(plot_violence_f4$names)

plot_violence_f4$namesfd <- factor(plot_violence_f4$namesf, levels = c("Believing Authorities \nthat Protesters \nwere Violent","Preemptive \nViolence \nPerception"))

#GROUP
plot_violence_f4$group <- plot_violence_f4$contrast

plot_violence_f4$group_f<-as.factor(plot_violence_f4$group)

plot_violence_f4$group_f<-dplyr::recode(plot_violence_f4$group_f, "2-1"="Generic Women \nvs. Men", "3-1"="Patriarchy-Compliant \nWomen vs. Men", "4-1"="Patriarchy-Defiant \nWomen vs. Men", "4-3"="Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women")

plot_violence_f4$group_f <- factor(plot_violence_f4$group_f, levels = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "Patriarchy-Defiant \nWomen vs. Men", "Patriarchy-Compliant \nWomen vs. Men", "Generic Women \nvs. Men"))


#PLOTTING
plot_v1v2_f4<-ggplot(plot_violence_f4, aes(x = estimate, y = namesfd, group=group_f, color=group_f))+ geom_point(aes(shape=group_f, color=group_f), size=3, position = position_dodge(width = 0.50))+ geom_errorbar(aes(xmin=conf.low, xmax=conf.high, group=group_f), size=1, width=0.3,position = position_dodge(width = 0.50)) + geom_vline(xintercept = 0,lty=2) + theme_bw()+theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5), axis.text.y = element_text(color = "grey20", size = 12, face = "bold"),legend.text=element_text(size=12, face = "bold")) + xlab("") + ylab("")


pl_f4 <- plot_v1v2_f4 + labs(color = " ", shape=" ") +scale_colour_grey(start = .15, end = .85)

p_f4 <- pl_f4 + theme(
  legend.position = "bottom", 
  legend.justification = "right", # Move the legend to the left
  legend.margin = margin(r = 7, unit = "mm") # Adjust this value to control the distance from the left side
)

p_f4

p_f4_c <- p_f4 +
  scale_shape_manual(values = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women" = 15, 
                                "Patriarchy-Defiant \nWomen vs. Men" = 17,  
                                "Patriarchy-Compliant \nWomen vs. Men" = 4,  
                                "Generic Women \nvs. Men" = 16)) 



# Reverse the order of the legend items
plot_v1v2_f4 <- p_f4_c + guides(color = guide_legend(reverse = TRUE),
                                     shape = guide_legend(reverse = TRUE))


plot_v1v2_f4



ggsave("GP_Figure4.tiff", plot_v1v2_f4, width = 9, height =6, dpi = 1000,
       units = "in")
```



#Figure 6 (Differences of Means for Immorality Perception)
```{r}
library(broom)
mhsd_df<- tidy(mhsd)
```

```{r}
mhsd_df_2 <- mhsd_df[c("contrast", "estimate","conf.low","conf.high", "adj.p.value")]

mhsd_df_2$dv <- "1"

mhsd_df_2_t <-mhsd_df_2[!(mhsd_df_2$contrast=="3-2" | mhsd_df_2$contrast=="4-2"),]

i1 <- mhsd_df_2_t$contrast=="4-3"
mhsd_df_2_t$estimate[i1] <- mhsd_df_2_t$estimate[i1] * (-1)
mhsd_df_2_t$conf.low[i1] <- mhsd_df_2_t$conf.low[i1] * (-1)
mhsd_df_2_t$conf.high[i1] <- mhsd_df_2_t$conf.high[i1] * (-1)
```


```{r}
plot_morality <- mhsd_df_2_t

```

```{r}
plot_morality$names <- plot_morality$dv

plot_morality$group <- plot_morality$contrast

plot_morality$group_f<-as.factor(plot_morality$group)

plot_morality$group_f<-dplyr::recode(plot_morality$group_f, "2-1"="Generic Women \nvs. Men", "3-1"="Patriarchy-Compliant \nWomen vs. Men", "4-1"="Patriarchy-Defiant \nWomen vs. Men", "4-3"="Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women")


plot_morality$group_fd <- factor(plot_morality$group_f, levels = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women","Patriarchy-Defiant \nWomen vs. Men","Patriarchy-Compliant \nWomen vs. Men", "Generic Women \nvs. Men"))


ploty1<-ggplot(plot_morality, aes(x = estimate, y = group_fd))+ geom_point()+ geom_errorbar(aes(xmin=conf.low, xmax=conf.high), size=1.15, width=0.15) + geom_vline(xintercept = 0,lty=2) + labs(title = "", subtitle = "") + theme_bw()+theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5)) + xlab("") + ylab("") 


plot11 <- ploty1 + theme(axis.text=element_text(size=15, face="bold"))
plot11


ggsave("GP_Figure6.tiff", plot11, width = 10, height = 5, dpi = 1000,
       units = "in")
```

#MEDIATION ANALYSES

#Create new dataframes 
#one that has only conditions 3 & 4 (Patriarchy-Compliant & Patriarchy-Defiant Women)
#one that has only conditions 1 & 3 (Men & Patriarchy-Compliant Women)
#one that has only conditions 1 & 4 (Men & Patriarchy-Defiant Women)
#one that has only conditions 1 & 2 (Men & Generic Women)


```{r}
library(tidyr)

gp_tm <- subset(gp_1, condition_f!=2)
gp_tm <- subset(gp_tm, condition_f!=4)
gp_tm$condition_f
gp_tm$cond_3 <- NA
gp_tm$cond_3 <- ifelse(gp_tm$condition_f == 3, 1, 0) # going from men (0) to traditional women (1)

gp_tn <- subset(gp_1, condition_f!=1)
gp_tn <- subset(gp_tn, condition_f!=2)
gp_tn$condition_f
gp_tn$cond_3 <- NA
gp_tn$cond_3 <- ifelse(gp_tn$condition_f == 3, 1, 0) # going from non-traditional (0) to traditional women (1)

gp_tn <- gp_tn %>% drop_na(toler_1_index)

gp_ntm <- subset(gp_1, condition_f!=2)
gp_ntm <- subset(gp_ntm, condition_f!=3)
gp_ntm$condition_f
gp_ntm$cond_3 <- NA
gp_ntm$cond_3 <- ifelse(gp_ntm$condition_f == 4, 1, 0) # going from men (0) to non-traditional women (1)

gp_ntm <- gp_ntm %>% drop_na(toler_1_index)


gp_gm <- subset(gp_1, condition_f!=3)
gp_gm <- subset(gp_gm, condition_f!=4)
gp_gm$condition_f
gp_gm$cond_3 <- NA
gp_gm$cond_3 <- ifelse(gp_gm$condition_f == 2, 1, 0) # going from men (0) to generic women (1)
```


#PERCEIVED VIOLENCE AS A MEDIATOR
#The below script can be used to replicate Figure 5 in the main body of the manuscript as well as the results in the "The Role of Gender in Countering Government Narratives" sub-section; Mediation analysis results are also provided in table form in SI-7 & SI-8, and the code for replicating these tables can be found in the "GP_SI_Analysis.Rmd"

# VIOLENCE 1 AS A MEDIATOR

```{r}
library("mediation")

set.seed(2014)
med.fit_v1_tn <- lm(violence_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)
out.fit_v1_tn <- lm(toler_1_index ~ violence_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)

med.out_v1_tn <- mediate(med.fit_v1_tn, out.fit_v1_tn, treat = "cond_3", mediator = "violence_1", boot = TRUE, sims = 1000)

summary(med.out_v1_tn)
#plot(med.out_v1_tn)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v1_tn <- med.out_v1_tn$d0
ACME_ci_lower_v1_tn <- med.out_v1_tn$d0.ci[1]
ACME_ci_upper_v1_tn <- med.out_v1_tn$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v1_tn <- med.out_v1_tn$z0
ADE_ci_lower_v1_tn <- med.out_v1_tn$z0.ci[1]
ADE_ci_upper_v1_tn <- med.out_v1_tn$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v1_tn <- med.out_v1_tn$tau.coef
TotalEffect_ci_lower_v1_tn <- med.out_v1_tn$tau.ci[1]
TotalEffect_ci_upper_v1_tn <- med.out_v1_tn$tau.ci[2]

summary_df_v1_tn <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v1_tn, ADE_estimate_v1_tn, TotalEffect_estimate_v1_tn),
  lwr = c(ACME_ci_lower_v1_tn, ADE_ci_lower_v1_tn, TotalEffect_ci_lower_v1_tn),
  upr = c(ACME_ci_upper_v1_tn, ADE_ci_upper_v1_tn, TotalEffect_ci_upper_v1_tn),
   group = "TN",
  dv = 1
)

```



```{r}
library("mediation")

set.seed(2014)
med.fit_v1_tm <- lm(violence_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)
out.fit_v1_tm <- lm(toler_1_index ~ violence_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)

med.out_v1_tm <- mediate(med.fit_v1_tm, out.fit_v1_tm, treat = "cond_3", mediator = "violence_1", boot = TRUE, sims = 1000)

summary(med.out_v1_tm)
#plot(med.out_v1_tm)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v1_tm <- med.out_v1_tm$d0
ACME_ci_lower_v1_tm <- med.out_v1_tm$d0.ci[1]
ACME_ci_upper_v1_tm <- med.out_v1_tm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v1_tm <- med.out_v1_tm$z0
ADE_ci_lower_v1_tm <- med.out_v1_tm$z0.ci[1]
ADE_ci_upper_v1_tm <- med.out_v1_tm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v1_tm <- med.out_v1_tm$tau.coef
TotalEffect_ci_lower_v1_tm <- med.out_v1_tm$tau.ci[1]
TotalEffect_ci_upper_v1_tm <- med.out_v1_tm$tau.ci[2]

summary_df_v1_tm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v1_tm, ADE_estimate_v1_tm, TotalEffect_estimate_v1_tm),
  lwr = c(ACME_ci_lower_v1_tm, ADE_ci_lower_v1_tm, TotalEffect_ci_lower_v1_tm),
  upr = c(ACME_ci_upper_v1_tm, ADE_ci_upper_v1_tm, TotalEffect_ci_upper_v1_tm),
   group = "TM",
  dv = 1
)

```



```{r}
library("mediation")

set.seed(2014)
med.fit_v1_gm <- lm(violence_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_gm)
out.fit_v1_gm <- lm(toler_1_index ~ violence_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_gm)

med.out_v1_gm <- mediate(med.fit_v1_gm, out.fit_v1_gm, treat = "cond_3", mediator = "violence_1", boot = TRUE, sims = 1000)

summary(med.out_v1_gm)
#plot(med.out_v1_gm)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v1_gm <- med.out_v1_gm$d0
ACME_ci_lower_v1_gm <- med.out_v1_gm$d0.ci[1]
ACME_ci_upper_v1_gm <- med.out_v1_gm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v1_gm <- med.out_v1_gm$z0
ADE_ci_lower_v1_gm <- med.out_v1_gm$z0.ci[1]
ADE_ci_upper_v1_gm <- med.out_v1_gm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v1_gm <- med.out_v1_gm$tau.coef
TotalEffect_ci_lower_v1_gm <- med.out_v1_gm$tau.ci[1]
TotalEffect_ci_upper_v1_gm <- med.out_v1_gm$tau.ci[2]

summary_df_v1_gm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v1_gm, ADE_estimate_v1_gm, TotalEffect_estimate_v1_gm),
  lwr = c(ACME_ci_lower_v1_gm, ADE_ci_lower_v1_gm, TotalEffect_ci_lower_v1_gm),
  upr = c(ACME_ci_upper_v1_gm, ADE_ci_upper_v1_gm, TotalEffect_ci_upper_v1_gm),
   group = "GM",
  dv = 1
)

```



```{r}
library("mediation")

set.seed(2014)
med.fit_v1_ntm <- lm(violence_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)
out.fit_v1_ntm <- lm(toler_1_index ~ violence_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)

med.out_v1_ntm <- mediate(med.fit_v1_ntm, out.fit_v1_ntm, treat = "cond_3", mediator = "violence_1", boot = TRUE, sims = 1000)

summary(med.out_v1_ntm)
#plot(med.out_v1_ntm)

```



```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v1_ntm <- med.out_v1_ntm$d0
ACME_ci_lower_v1_ntm <- med.out_v1_ntm$d0.ci[1]
ACME_ci_upper_v1_ntm <- med.out_v1_ntm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v1_ntm <- med.out_v1_ntm$z0
ADE_ci_lower_v1_ntm <- med.out_v1_ntm$z0.ci[1]
ADE_ci_upper_v1_ntm <- med.out_v1_ntm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v1_ntm <- med.out_v1_ntm$tau.coef
TotalEffect_ci_lower_v1_ntm <- med.out_v1_ntm$tau.ci[1]
TotalEffect_ci_upper_v1_ntm <- med.out_v1_ntm$tau.ci[2]

summary_df_v1_ntm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v1_ntm, ADE_estimate_v1_ntm, TotalEffect_estimate_v1_ntm),
  lwr = c(ACME_ci_lower_v1_ntm, ADE_ci_lower_v1_ntm, TotalEffect_ci_lower_v1_ntm),
  upr = c(ACME_ci_upper_v1_ntm, ADE_ci_upper_v1_ntm, TotalEffect_ci_upper_v1_ntm),
   group = "NTM",
  dv = 1
)

```


# VIOLENCE 2 AS A MEDIATOR

```{r}
library("mediation")

set.seed(2014)
med.fit_v2_tn <- lm(violence_2 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox +uni, data = gp_tn)
out.fit_v2_tn <- lm(toler_2_index ~ violence_2 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox +uni, data = gp_tn)

med.out_v2_tn <- mediate(med.fit_v2_tn, out.fit_v2_tn, treat = "cond_3", mediator = "violence_2", boot = TRUE, sims = 1000)

summary(med.out_v2_tn)
#plot(med.out_v2_tn)

```


```{r}
sens.out_v2_tn <- medsens(med.out_v2_tn, rho.by = 0.1, effect.type = "indirect", sims = 1000)
summary(sens.out_v2_tn)

```


```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v2_tn <- med.out_v2_tn$d0
ACME_ci_lower_v2_tn <- med.out_v2_tn$d0.ci[1]
ACME_ci_upper_v2_tn <- med.out_v2_tn$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v2_tn <- med.out_v2_tn$z0
ADE_ci_lower_v2_tn <- med.out_v2_tn$z0.ci[1]
ADE_ci_upper_v2_tn <- med.out_v2_tn$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v2_tn <- med.out_v2_tn$tau.coef
TotalEffect_ci_lower_v2_tn <- med.out_v2_tn$tau.ci[1]
TotalEffect_ci_upper_v2_tn <- med.out_v2_tn$tau.ci[2]

summary_df_v2_tn <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v2_tn, ADE_estimate_v2_tn, TotalEffect_estimate_v2_tn),
  lwr = c(ACME_ci_lower_v2_tn, ADE_ci_lower_v2_tn, TotalEffect_ci_lower_v2_tn),
  upr = c(ACME_ci_upper_v2_tn, ADE_ci_upper_v2_tn, TotalEffect_ci_upper_v2_tn),
   group = "TN",
  dv = 2
)

```



```{r}
library("mediation")

set.seed(2014)
med.fit_v2_tm <- lm(violence_2 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox +uni, data = gp_tm)
out.fit_v2_tm <- lm(toler_2_index ~ violence_2 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox +uni, data = gp_tm)

med.out_v2_tm <- mediate(med.fit_v2_tm, out.fit_v2_tm, treat = "cond_3", mediator = "violence_2", boot = TRUE, sims = 1000)

summary(med.out_v2_tm)
#plot(med.out_v2_tm)

```



```{r}
sens.out_v2_tm <- medsens(med.out_v2_tm, rho.by = 0.1, effect.type = "indirect", sims = 1000)
summary(sens.out_v2_tm)


```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_v2_tm <- med.out_v2_tm$d0
ACME_ci_lower_v2_tm <- med.out_v2_tm$d0.ci[1]
ACME_ci_upper_v2_tm <- med.out_v2_tm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_v2_tm <- med.out_v2_tm$z0
ADE_ci_lower_v2_tm <- med.out_v2_tm$z0.ci[1]
ADE_ci_upper_v2_tm <- med.out_v2_tm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_v2_tm <- med.out_v2_tm$tau.coef
TotalEffect_ci_lower_v2_tm <- med.out_v2_tm$tau.ci[1]
TotalEffect_ci_upper_v2_tm <- med.out_v2_tm$tau.ci[2]

summary_df_v2_tm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_v2_tm, ADE_estimate_v2_tm, TotalEffect_estimate_v2_tm),
  lwr = c(ACME_ci_lower_v2_tm, ADE_ci_lower_v2_tm, TotalEffect_ci_lower_v2_tm),
  upr = c(ACME_ci_upper_v2_tm, ADE_ci_upper_v2_tm, TotalEffect_ci_upper_v2_tm),
   group = "TM",
  dv = 2
)

```


```{r}
summary_df_combined_v1v2_mediation <- rbind(summary_df_v1_tn,summary_df_v1_tm,summary_df_v1_gm,summary_df_v1_ntm,summary_df_v2_tm,summary_df_v2_tn)

#write.csv(summary_df_combined_v1v2_mediation, file = "v1v2_mediation_results_2023.csv", row.names = FALSE)
```

```{r}
Plotting_violence_mediation <- summary_df_combined_v1v2_mediation

```

#Figure 5: The Mediating Effect of Perceived Violence on Citizens’ Tolerance of Repression
```{r}
library(ggplot2)
Plotting_violence_mediation$names <- Plotting_violence_mediation$dv

Plotting_violence_mediation$names<-dplyr::recode(Plotting_violence_mediation$names, 
                          "1"="A) Preemptive Violence Perception & \nRepression Toleration",
                          "2"="B) Reactive Violence Perception & \nRepression Toleration")

Plotting_violence_mediation$Effect<-dplyr::recode(Plotting_violence_mediation$Effect, 
                          "Total Effect"="Total \nEffect")

Plotting_violence_mediation$namesf<-as.factor(Plotting_violence_mediation$names)

Plotting_violence_mediation$namesfd <- factor(Plotting_violence_mediation$namesf, levels = c("A) Preemptive Violence Perception & \nRepression Toleration","B) Reactive Violence Perception & \nRepression Toleration"))

Plotting_violence_mediation$effect_name <- factor(Plotting_violence_mediation$Effect, levels = c("Total \nEffect","ADE","ACME"))

Plotting_violence_mediation$group <-dplyr::recode(Plotting_violence_mediation$group, "TN"="Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "TM"="Patriarchy-Compliant \nWomen vs. Men","NTM"="Patriarchy-Defiant \nWomen vs. Men","GM"="Generic Women \nvs. Men")

Plotting_violence_mediation$group_f<-as.factor(Plotting_violence_mediation$group)

Plotting_violence_mediation$group_f <- factor(Plotting_violence_mediation$group_f, levels = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "Patriarchy-Defiant \nWomen vs. Men", "Patriarchy-Compliant \nWomen vs. Men", "Generic Women \nvs. Men"))


plot_mediation_v_t1t2<-ggplot(Plotting_violence_mediation, aes(x = Estimate, y = effect_name, group=group_f, color=group_f))+ geom_point(aes(shape=group_f, color=group_f), size=2, position = position_dodge(width = 0.50))+ geom_errorbar(aes(xmin=lwr, xmax=upr, group=group_f), size=1, width=0.3,position = position_dodge(width = 0.50)) + geom_vline(xintercept = 0,lty=2) + labs(title = "") + theme_bw()+theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5), axis.text.y = element_text(color = "grey20", size = 14, face = "bold"),legend.text=element_text(color = "black", size=12, face = "bold")) + xlab("") + ylab("") + facet_wrap(~namesfd,ncol=2) + theme(strip.text = element_text(color = "black",size=14, face = "bold"))


plot_mediation_v_t1t2 <- plot_mediation_v_t1t2 + labs(color = " ", shape=" ") + scale_colour_grey(start = .15, end = .85)

plot_mediation_v_t1t2 <- plot_mediation_v_t1t2 + theme(
  legend.position = "bottom", 
  legend.justification = "right", # Move the legend to the left
  legend.margin = margin(r = 30, unit = "mm"))


plot_mediation_v_t1t2_2c <- plot_mediation_v_t1t2 + theme(strip.background = element_rect(fill = gray(0.92), color = "black"), strip.placement = "free")

plot_mediation_v_t1t2_2c <- plot_mediation_v_t1t2_2c +
  scale_shape_manual(values = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women" = 15, 
                                "Patriarchy-Defiant \nWomen vs. Men" = 17,  
                                "Patriarchy-Compliant \nWomen vs. Men" = 4,  
                                "Generic Women \nvs. Men" = 16)) 



# Reverse the order of the legend items
plot_mediation_v_t1t2_f5 <- plot_mediation_v_t1t2_2c + guides(color = guide_legend(reverse = TRUE),
                                     shape = guide_legend(reverse = TRUE))


plot_mediation_v_t1t2_f5


ggsave("GP_Figure5.tiff", plot_mediation_v_t1t2_f5, width = 11, height = 6, dpi = 1000,
       units = "in")

```


# PERCEIVED MORALITY AS A MEDIATOR

#The below script can be used to replicate Figure 7 in the main body of the manuscript as well as the results in the "The Role of Perceived Morality" sub-section; Mediation analysis results are also provided in table form in SI-7 & SI-8, and the code for replicating these tables can be found in the "GP_SI_Analysis.Rmd"

#Tolerance 1
```{r}
library("mediation")

set.seed(2014)
med.fit_m_t1_tn <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)
out.fit_m_t1_tn <- lm(toler_1_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)

med.out_m_t1_tn <- mediate(med.fit_m_t1_tn, out.fit_m_t1_tn, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t1_tn)
#plot(med.out_m_t1_tn)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate <- med.out_m_t1_tn$d0
ACME_ci_lower <- med.out_m_t1_tn$d0.ci[1]
ACME_ci_upper <- med.out_m_t1_tn$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate <- med.out_m_t1_tn$z0
ADE_ci_lower <- med.out_m_t1_tn$z0.ci[1]
ADE_ci_upper <- med.out_m_t1_tn$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate <- med.out_m_t1_tn$tau.coef
TotalEffect_ci_lower <- med.out_m_t1_tn$tau.ci[1]
TotalEffect_ci_upper <- med.out_m_t1_tn$tau.ci[2]

summary_df_m_t1_tn <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate, ADE_estimate, TotalEffect_estimate),
  lwr = c(ACME_ci_lower, ADE_ci_lower, TotalEffect_ci_lower),
  upr = c(ACME_ci_upper, ADE_ci_upper, TotalEffect_ci_upper),
   group = "TNT",
  dv = 1
)

```




```{r}
library("mediation")

set.seed(2014)
med.fit_m_t1_tm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)
out.fit_m_t1_tm <- lm(toler_1_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)

med.out_m_t1_tm <- mediate(med.fit_m_t1_tm, out.fit_m_t1_tm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t1_tm)
#plot(med.out_m_t1_tm)

```




```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_tm1 <- med.out_m_t1_tm$d0
ACME_ci_lower_tm1 <- med.out_m_t1_tm$d0.ci[1]
ACME_ci_upper_tm1 <- med.out_m_t1_tm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_tm1 <- med.out_m_t1_tm$z0
ADE_ci_lower_tm1 <- med.out_m_t1_tm$z0.ci[1]
ADE_ci_upper_tm1 <- med.out_m_t1_tm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_tm1 <- med.out_m_t1_tm$tau.coef
TotalEffect_ci_lower_tm1 <- med.out_m_t1_tm$tau.ci[1]
TotalEffect_ci_upper_tm1 <- med.out_m_t1_tm$tau.ci[2]

summary_df_m_t1_tm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_tm1, ADE_estimate_tm1, TotalEffect_estimate_tm1),
  lwr = c(ACME_ci_lower_tm1, ADE_ci_lower_tm1, TotalEffect_ci_lower_tm1),
  upr = c(ACME_ci_upper_tm1, ADE_ci_upper_tm1, TotalEffect_ci_upper_tm1),
   group = "TM",
  dv = 1
)

```


```{r}
library("mediation")

set.seed(2014)
med.fit_m_t1_gm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_gm)
out.fit_m_t1_gm <- lm(toler_1_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_gm)

med.out_m_t1_gm <- mediate(med.fit_m_t1_gm, out.fit_m_t1_gm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t1_gm)
#plot(med.out_m_t1_gm)

```


```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_gm1 <- med.out_m_t1_gm$d0
ACME_ci_lower_gm1 <- med.out_m_t1_gm$d0.ci[1]
ACME_ci_upper_gm1 <- med.out_m_t1_gm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_gm1 <- med.out_m_t1_gm$z0
ADE_ci_lower_gm1 <- med.out_m_t1_gm$z0.ci[1]
ADE_ci_upper_gm1 <- med.out_m_t1_gm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_gm1 <- med.out_m_t1_gm$tau.coef
TotalEffect_ci_lower_gm1 <- med.out_m_t1_gm$tau.ci[1]
TotalEffect_ci_upper_gm1 <- med.out_m_t1_gm$tau.ci[2]

summary_df_m_t1_gm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_gm1, ADE_estimate_gm1, TotalEffect_estimate_gm1),
  lwr = c(ACME_ci_lower_gm1, ADE_ci_lower_gm1, TotalEffect_ci_lower_gm1),
  upr = c(ACME_ci_upper_gm1, ADE_ci_upper_gm1, TotalEffect_ci_upper_gm1),
   group = "GM",
  dv = 1
)

```


```{r}
library("mediation")

set.seed(2014)
med.fit_m_t1_ntm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)
out.fit_m_t1_ntm <- lm(toler_1_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)

med.out_m_t1_ntm <- mediate(med.fit_m_t1_ntm, out.fit_m_t1_ntm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t1_ntm)
#plot(med.out_m_t1_ntm)

```


```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_ntm1 <- med.out_m_t1_ntm$d0
ACME_ci_lower_ntm1 <- med.out_m_t1_ntm$d0.ci[1]
ACME_ci_upper_ntm1 <- med.out_m_t1_ntm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_ntm1 <- med.out_m_t1_ntm$z0
ADE_ci_lower_ntm1 <- med.out_m_t1_ntm$z0.ci[1]
ADE_ci_upper_ntm1 <- med.out_m_t1_ntm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_ntm1 <- med.out_m_t1_ntm$tau.coef
TotalEffect_ci_lower_ntm1 <- med.out_m_t1_ntm$tau.ci[1]
TotalEffect_ci_upper_ntm1 <- med.out_m_t1_ntm$tau.ci[2]

summary_df_m_t1_ntm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_ntm1, ADE_estimate_ntm1, TotalEffect_estimate_ntm1),
  lwr = c(ACME_ci_lower_ntm1, ADE_ci_lower_ntm1, TotalEffect_ci_lower_ntm1),
  upr = c(ACME_ci_upper_ntm1, ADE_ci_upper_ntm1, TotalEffect_ci_upper_ntm1),
   group = "MNT",
  dv = 1
)

```


#Tolerance 2

```{r}
library("mediation")

set.seed(2014)
med.fit_m_t2_tn <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)
out.fit_m_t2_tn <- lm(toler_2_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)

med.out_m_t2_tn <- mediate(med.fit_m_t2_tn, out.fit_m_t2_tn, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t2_tn)
#plot(med.out_m_t2_tn)

```


```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_tn2 <- med.out_m_t2_tn$d0
ACME_ci_lower_tn2 <- med.out_m_t2_tn$d0.ci[1]
ACME_ci_upper_tn2 <- med.out_m_t2_tn$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_tn2 <- med.out_m_t2_tn$z0
ADE_ci_lower_tn2 <- med.out_m_t2_tn$z0.ci[1]
ADE_ci_upper_tn2 <- med.out_m_t2_tn$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_tn2 <- med.out_m_t2_tn$tau.coef
TotalEffect_ci_lower_tn2 <- med.out_m_t2_tn$tau.ci[1]
TotalEffect_ci_upper_tn2 <- med.out_m_t2_tn$tau.ci[2]

summary_df_m_t2_tn <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_tn2, ADE_estimate_tn2, TotalEffect_estimate_tn2),
  lwr = c(ACME_ci_lower_tn2, ADE_ci_lower_tn2, TotalEffect_ci_lower_tn2),
  upr = c(ACME_ci_upper_tn2, ADE_ci_upper_tn2, TotalEffect_ci_upper_tn2),
   group = "TNT",
  dv = 2
)

```


```{r}
library("mediation")

set.seed(2014)
med.fit_m_t2_tm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)
out.fit_m_t2_tm <- lm(toler_2_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)

med.out_m_t2_tm <- mediate(med.fit_m_t2_tm, out.fit_m_t2_tm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t2_tm)
#plot(med.out_m_t2_tm)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_tm2 <- med.out_m_t2_tm$d0
ACME_ci_lower_tm2 <- med.out_m_t2_tm$d0.ci[1]
ACME_ci_upper_tm2 <- med.out_m_t2_tm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_tm2 <- med.out_m_t2_tm$z0
ADE_ci_lower_tm2 <- med.out_m_t2_tm$z0.ci[1]
ADE_ci_upper_tm2 <- med.out_m_t2_tm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_tm2 <- med.out_m_t2_tm$tau.coef
TotalEffect_ci_lower_tm2 <- med.out_m_t2_tm$tau.ci[1]
TotalEffect_ci_upper_tm2 <- med.out_m_t2_tm$tau.ci[2]

summary_df_m_t2_tm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_tm2, ADE_estimate_tm2, TotalEffect_estimate_tm2),
  lwr = c(ACME_ci_lower_tm2, ADE_ci_lower_tm2, TotalEffect_ci_lower_tm2),
  upr = c(ACME_ci_upper_tm2, ADE_ci_upper_tm2, TotalEffect_ci_upper_tm2),
   group = "TM",
  dv = 2
)

```




```{r}
library("mediation")

set.seed(2014)
med.fit_m_t2_ntm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)
out.fit_m_t2_ntm <- lm(toler_2_index ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)

med.out_m_t2_ntm <- mediate(med.fit_m_t2_ntm, out.fit_m_t2_ntm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_t2_ntm)
#plot(med.out_m_t2_ntm)

```


```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_ntm2 <- med.out_m_t2_ntm$d0
ACME_ci_lower_ntm2 <- med.out_m_t2_ntm$d0.ci[1]
ACME_ci_upper_ntm2 <- med.out_m_t2_ntm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_ntm2 <- med.out_m_t2_ntm$z0
ADE_ci_lower_ntm2 <- med.out_m_t2_ntm$z0.ci[1]
ADE_ci_upper_ntm2 <- med.out_m_t2_ntm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_ntm2 <- med.out_m_t2_ntm$tau.coef
TotalEffect_ci_lower_ntm2 <- med.out_m_t2_ntm$tau.ci[1]
TotalEffect_ci_upper_ntm2 <- med.out_m_t2_ntm$tau.ci[2]

summary_df_m_t2_ntm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_ntm2, ADE_estimate_ntm2, TotalEffect_estimate_ntm2),
  lwr = c(ACME_ci_lower_ntm2, ADE_ci_lower_ntm2, TotalEffect_ci_lower_ntm2),
  upr = c(ACME_ci_upper_ntm2, ADE_ci_upper_ntm2, TotalEffect_ci_upper_ntm2),
   group = "MNT",
  dv = 2
)

```


# Believing Authorities

```{r}
library("mediation")

set.seed(2014)
med.fit_m_v_tn <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)
out.fit_m_v_tn <- lm(violence_2 ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tn)

med.out_m_v_tn <- mediate(med.fit_m_v_tn, out.fit_m_v_tn, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_v_tn)
#plot(med.out_m_v_tn)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_tn3 <- med.out_m_v_tn$d0
ACME_ci_lower_tn3 <- med.out_m_v_tn$d0.ci[1]
ACME_ci_upper_tn3 <- med.out_m_v_tn$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_tn3 <- med.out_m_v_tn$z0
ADE_ci_lower_tn3 <- med.out_m_v_tn$z0.ci[1]
ADE_ci_upper_tn3 <- med.out_m_v_tn$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_tn3 <- med.out_m_v_tn$tau.coef
TotalEffect_ci_lower_tn3 <- med.out_m_v_tn$tau.ci[1]
TotalEffect_ci_upper_tn3 <- med.out_m_v_tn$tau.ci[2]

summary_df_m_v_tn <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_tn3, ADE_estimate_tn3, TotalEffect_estimate_tn3),
  lwr = c(ACME_ci_lower_tn3, ADE_ci_lower_tn3, TotalEffect_ci_lower_tn3),
  upr = c(ACME_ci_upper_tn3, ADE_ci_upper_tn3, TotalEffect_ci_upper_tn3),
   group = "TNT",
  dv = 3
)


```




```{r}
library("mediation")

set.seed(2014)
med.fit_m_v_tm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)
out.fit_m_v_tm <- lm(violence_2 ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_tm)

med.out_m_v_tm <- mediate(med.fit_m_v_tm, out.fit_m_v_tm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_v_tm)
#plot(med.out_m_v_tm)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_tm3 <- med.out_m_v_tm$d0
ACME_ci_lower_tm3 <- med.out_m_v_tm$d0.ci[1]
ACME_ci_upper_tm3 <- med.out_m_v_tm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_tm3 <- med.out_m_v_tm$z0
ADE_ci_lower_tm3 <- med.out_m_v_tm$z0.ci[1]
ADE_ci_upper_tm3 <- med.out_m_v_tm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_tm3 <- med.out_m_v_tm$tau.coef
TotalEffect_ci_lower_tm3 <- med.out_m_v_tm$tau.ci[1]
TotalEffect_ci_upper_tm3 <- med.out_m_v_tm$tau.ci[2]

summary_df_m_v_tm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_tm3, ADE_estimate_tm3, TotalEffect_estimate_tm3),
  lwr = c(ACME_ci_lower_tm3, ADE_ci_lower_tm3, TotalEffect_ci_lower_tm3),
  upr = c(ACME_ci_upper_tm3, ADE_ci_upper_tm3, TotalEffect_ci_upper_tm3),
   group = "TM",
  dv = 3
)


```



```{r}
library("mediation")

set.seed(2014)
med.fit_m_v_ntm <- lm(mor_1 ~ cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)
out.fit_m_v_ntm <- lm(violence_2 ~ mor_1 + cond_3 + age_b_st + sex_b + income_b_st + putin_support_st + orthodox + uni, data = gp_ntm)

med.out_m_v_ntm <- mediate(med.fit_m_v_ntm, out.fit_m_v_ntm, treat = "cond_3", mediator = "mor_1", boot = TRUE, sims = 1000)

summary(med.out_m_v_ntm)
#plot(med.out_m_v_ntm)

```

```{r}
# Extract the estimates and 95% Confidence Intervals for ACME
ACME_estimate_ntm3 <- med.out_m_v_ntm$d0
ACME_ci_lower_ntm3 <- med.out_m_v_ntm$d0.ci[1]
ACME_ci_upper_ntm3 <- med.out_m_v_ntm$d0.ci[2]

# Extract the estimates and 95% Confidence Intervals for ADE
ADE_estimate_ntm3 <- med.out_m_v_ntm$z0
ADE_ci_lower_ntm3 <- med.out_m_v_ntm$z0.ci[1]
ADE_ci_upper_ntm3 <- med.out_m_v_ntm$z0.ci[2]

# Extract the estimates and 95% Confidence Intervals for Total Effect
TotalEffect_estimate_ntm3 <- med.out_m_v_ntm$tau.coef
TotalEffect_ci_lower_ntm3 <- med.out_m_v_ntm$tau.ci[1]
TotalEffect_ci_upper_ntm3 <- med.out_m_v_ntm$tau.ci[2]

summary_df_m_v_ntm <- data.frame(
  Effect = c("ACME", "ADE", "Total Effect"),
  Estimate = c(ACME_estimate_ntm3, ADE_estimate_ntm3, TotalEffect_estimate_ntm3),
  lwr = c(ACME_ci_lower_ntm3, ADE_ci_lower_ntm3, TotalEffect_ci_lower_ntm3),
  upr = c(ACME_ci_upper_ntm3, ADE_ci_upper_ntm3, TotalEffect_ci_upper_ntm3),
   group = "MNT",
  dv = 3
)


```


```{r}
summary_df_combined <- rbind(summary_df_m_t1_tm, summary_df_m_t2_tm, summary_df_m_v_tm, 
                             summary_df_m_t1_tn, summary_df_m_t2_tn, summary_df_m_v_tn, 
                             summary_df_m_t1_ntm, summary_df_m_t2_ntm, summary_df_m_v_ntm)

#write.csv(summary_df_combined, file = "morality_plot_results_2023.csv", row.names = FALSE)
```


# Figure 7 (Morality as a Mediator)
```{r}
library(readxl)

Plotting_morality_mediation <- summary_df_combined
```


```{r}
library(ggplot2)
Plotting_morality_mediation$names <- Plotting_morality_mediation$dv

Plotting_morality_mediation$names<-dplyr::recode(Plotting_morality_mediation$names, 
                          "1"="A) Preemptive \nRepression Toleration",
                          "2"="B) Reactive \nRepression Toleration",
                          "3"="C) Believing Authorities' Claim \nabout Protest Violence")


Plotting_morality_mediation$namesf<-as.factor(Plotting_morality_mediation$names)

Plotting_morality_mediation$namesfd <- factor(Plotting_morality_mediation$namesf, levels = c("A) Preemptive \nRepression Toleration","B) Reactive \nRepression Toleration","C) Believing Authorities' Claim \nabout Protest Violence"))


Plotting_morality_mediation$Effect<-dplyr::recode(Plotting_morality_mediation$Effect, "Total Effect"="Total \nEffect")

Plotting_morality_mediation$effect_name <- factor(Plotting_morality_mediation$Effect, levels = c("Total \nEffect","ADE","ACME"))



Plotting_morality_mediation$group<-dplyr::recode(Plotting_morality_mediation$group, "TNT"="Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "TM"="Patriarchy-Compliant \nWomen vs. Men","MNT"="Patriarchy-Defiant \nWomen vs. Men")


Plotting_morality_mediation$group_f<-as.factor(Plotting_morality_mediation$group)

Plotting_morality_mediation$group_f <- factor(Plotting_morality_mediation$group_f, levels = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women", "Patriarchy-Defiant \nWomen vs. Men", "Patriarchy-Compliant \nWomen vs. Men"))


plot_mediation_m<-ggplot(Plotting_morality_mediation, aes(x = Estimate, y = effect_name, group=group_f, color=group_f)) + geom_point(aes(shape=group_f, color=group_f), size=2, position = position_dodge(width = 0.50))+ geom_errorbar(aes(xmin=lwr, xmax=upr, group=group_f), size=1, width=0.3,position = position_dodge(width = 0.50)) + geom_vline(xintercept = 0,lty=2) + labs(title = "") + theme_bw()+theme(plot.title = element_text(hjust = 0.5), plot.subtitle = element_text(hjust = 0.5), axis.text.y = element_text(color = "black", size = 15, face = "bold"),legend.text=element_text(color = "black", size=15, face = "bold")) + xlab("") + ylab("") + facet_wrap(~namesfd,ncol=3) + theme(strip.text = element_text(color = "black",size=16, face = "bold"))


plot_mediation_m <- plot_mediation_m + labs(color = " ", shape=" ") + scale_colour_grey(start = .10, end = .70)
plot_mediation_m <- plot_mediation_m + theme(
  legend.position = "bottom", 
  legend.justification = "right", # Move the legend to the left
  legend.margin = margin(r = 55, unit = "mm"))

plot_mediation_m_c <- plot_mediation_m + theme(strip.background = element_rect(fill = gray(0.92), color = "black"), strip.placement = "free")


plot_mediation_m_c_2c <- plot_mediation_m_c +
  scale_shape_manual(values = c("Patriarchy-Compliant Women \nvs. Patriarchy-Defiant Women" = 15, 
                                "Patriarchy-Defiant \nWomen vs. Men" = 17,  
                                "Patriarchy-Compliant \nWomen vs. Men" = 4)) 



# Reverse the order of the legend items
plot_mediation_m_c_2c_f7 <- plot_mediation_m_c_2c + guides(color = guide_legend(reverse = TRUE),
                                     shape = guide_legend(reverse = TRUE))


plot_mediation_m_c_2c_f7



ggsave("GP_Figure7.tiff", plot_mediation_m_c_2c_f7, width = 14, height = 6, dpi = 1000,
       units = "in")

```



#Subset Analysis by Respondent's Gender 

```{r}
#More Analyses on Treatment Heterogeneity in Appendix; This code can be used to replicate results listed in the third paragraph of the "Protesters’ Gender Influences Citizens’ Tolerance of Repression" subsection in the Results section
 

g_t1_tn_m <- aov(toler_1_index ~ condition_f, data = subset(gp_tn, sex_b==1))
TukeyHSD(g_t1_tn_m)

g_t1_tn_w <- aov(toler_1_index ~ condition_f, data = subset(gp_tn, sex_b==0))
TukeyHSD(g_t1_tn_w)

```

# Petition Signing

```{r}
gp_1$petition_sign_f <- as.factor(gp_1$petition_sign)
gp_1$petition_sign_f <- relevel(gp_1$petition_sign_f, ref="no signature")
gp_1 <- within(gp_1, condition_f <- relevel(condition_f, ref = 3))
```

```{r}
library(nnet)
test <- multinom(petition_sign_f ~ condition_f, data = gp_1)
summary(test)
```

```{r}
library(texreg)
texreg(l = list(test),stars = c(0.001, 0.01, 0.05))

```


```{r}
coefs <- coef(test)
coefs
```

```{r}
exp(coef(test))
```

```{r}
(exp(coefs)-1)*100
```

```{r}
z <- summary(test)$coefficients/summary(test)$standard.errors
z
p <- (1 - pnorm(abs(z), 0, 1)) * 2
p
```

#Cronbach's alphas 

```{r}
library(ltm)

toler_1_index_ca <- data.frame(gp_1$tolerance_1, gp_1$tolerance_2, gp_1$tolerance_3, gp_1$tolerance_4)
cronbach.alpha(toler_1_index_ca, na.rm = T) #0.873

toler_2_index_ca <- data.frame(gp_1$t_2_1, gp_1$t_2_2, gp_1$t_2_3, gp_1$t_2_4, gp_1$t_2_5)
cronbach.alpha(toler_2_index_ca, na.rm = T) #0.88

```




